clear all;
close all;
fclose all;

%% Register all processed HCR images

path1='C:\Users\Tom�s\Documents\1. Regehr Lab\Data\HCR\Registration and analysis\0. raw single Z\';
% get the folder contents
d2 = dir(path1); %using this directory to extract slice names is useful because it only contains bacthes 2 and 3 which are the ones that were processed in Matlab.
% remove all files (isdir property is 0)
dfolders = d2([d2(:).isdir]==1);
% remove '.' and '..'
dfolders = dfolders(~ismember({dfolders(:).name},{'.','..'}));


path2='C:\Users\Tom�s\Documents\1. Regehr Lab\Data\HCR\Registration and analysis\1.Processed not registered\'; %processed channels are here use find string or change file names to all same format.
cppath='C:\Users\Tom�s\Documents\1. Regehr Lab\Data\HCR\Registration and analysis\2.registered after downscaling to half res\';

%load control points
load([cppath 'controlpointsbatch2.mat'])
load([cppath 'controlpointsbatch3.mat'])

%load fixed reference
fixed=imread([path2 's4b1' '\'  'REFERENCEHalfRes.tif']); %matlab format out for tif is .tiff

    f=waitbar(0,'starting')
tic
for i=1:numel(dfolders)
   
     fracprogress=(i-1)/numel(dfolders);
      waitbar(fracprogress, f, {['Registering HCR images'];['Remaining time: ' num2str((1-fracprogress)*(toc/fracprogress)/60) ' minutes']} );
            
    
    %use eval and stringcat to know correct cp variable names.
    currentSlice=dfolders(i).name;
    fixedPointsVarName=['fixedPoints' currentSlice(2:end) '_4b1'];
    movingPointsVarName=['movingPoints' currentSlice(2:end) '_4b1'];
    currentFixedPoints= eval(fixedPointsVarName);
     currentMovingPoints=eval(movingPointsVarName);

    tform=fitgeotrans(currentMovingPoints,currentFixedPoints,'polynomial',4);
    
    %load 3 channels and background
     background=imread([path2 dfolders(i).name '\'  dfolders(i).name 'backgroundHalfRes.tiff']); %matlab format out for tif is .tiff
     Nxph1=imread([path2 dfolders(i).name '\'  dfolders(i).name 'Nxph1HalfRes.tiff']); %matlab format out for tif is .tiff
     Slc6a5=imread([path2 dfolders(i).name '\'  dfolders(i).name 'Slc6a5HalfRes.tiff']); %matlab format out for tif is .tiff
     Aldh1a3=imread([path2 dfolders(i).name '\'  dfolders(i).name 'Aldh1a3HalfRes.tiff']); %matlab format out for tif is .tiff
     
     %register images to fixed reference (s4b1)
     backgroundRegistered=imwarp(background,tform,'OutputView',imref2d(size(fixed)));
     Nxph1Registered=imwarp(Nxph1,tform,'OutputView',imref2d(size(fixed)));
     Slc6a5Registered=imwarp(Slc6a5,tform,'OutputView',imref2d(size(fixed)));
     Aldh1a3Registered=imwarp(Aldh1a3,tform,'OutputView',imref2d(size(fixed)));
     
     % Save all registered images
     imwrite(backgroundRegistered,[cppath currentSlice '\' 'registered_' currentSlice '_background.tiff']);
     imwrite(Nxph1Registered,[cppath currentSlice '\' 'registered_' currentSlice '_Nxph1.tiff']);
     imwrite(Slc6a5Registered,[cppath currentSlice '\' 'registered_' currentSlice '_Slc6a5.tiff']);
     imwrite(Aldh1a3Registered,[cppath currentSlice '\' 'registered_' currentSlice '_Aldh1a3.tiff']);
     
     figure;
     imshowpair(fixed,backgroundRegistered)
     title(currentSlice)
     
end
     waitbar(1, f,'Done!')
     toc

       


% % figure
% % imshowpair(fixed*10,movingregistered54*10)
% 
% 
% % imwrite(registeredSlice5Aldhla3,'regslice5.tiff')
% % 
% currentSlice='s2b3'
% tform=fitgeotrans(movingPoints2b3_4b1,fixedPoints2b3_4b1,'polynomial',4);
% registered=imwarp(Nxph1,tform,'OutputView',imref2d(size(fixed)));
% 
%      imwrite(registered,[cppath currentSlice '\' 'registered_' currentSlice '_Nxph1.tiff']);


% % registered34_slice3_Aldh1a3=imwarp(slice3Aldhla3,tform,'OutputView',imref2d(size(fixed)));
% 
% 
% movingregistered=imwarp(moving,tform,'OutputView',imref2d(size(fixed)));
% imwrite(movingregistered,[path1 'registered_4b2_4b1_Slc6a5.tiff'])
% % imwrite(registered54_slice5_Aldh1a3,'registered54_slice5_Aldh1a3.tiff')


